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PREFACE 


The  ever-increasing  demands  for  greater  measurement  accuracy  require  the  careful 
investigation  of  effects  previously  treated  heuristicaUy.  What  would  have  been  considered 
hairsplitting  yesterday  is  expected  today  and  may  be  taken  for  granted  tomorrow.  An 
example  is  the  delay  calibration  problem,  which  is  the  subject  of  this  little  monograph. 

For  a  full  appreciation  of  the  theory  presented,  the  reader  should  be  familiar  with 
probability  and  statistics  at  least  on  the  level  of  the  book  by  Freund  [1] .  The  practical 
application  of  the  statistical  methods  described  requires  the  extensive  use  of  computer 
programs.  The  book  by  Afifi  and  Azen  [2]  provides  a  valuable  exposition  on  the  selec¬ 
tion  and  interpretation  of  “packaged”  statistical  programs. 

1  have  attempted  to  describe  the  measurement  methods  in  a  general  way  so  that 
many  variations  can  be  included  in  the  same  general  framework.  The  mathematical 
models  are  thought  to  be  the  simplest  ones  adequate  to  describe  the  physical  reality. 
Despite  this  simplicity,  the  models  yield  a  surprising  richness  of  insight  and  mathematical 
content. 

The  measurement  methods,  the  statistical  techniques,  and  the  theoretical  results  have 
been  thoroughly  tested.  In  this  regard,  I  gratefully  acknowledge  the  assistance  of  Hans 
Ruhr,  Douglas  Wahrenburger,  and  Dean  Watkeys  of  Locus,  Inc.  Their  contributions 
include  the  assembling  of  instrumentation,  taking  of  data,  writing  of  programs,  and  the 
performance  of  Monte  Carlo  studies  which  helped  verify  the  theory. 
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METHODS,  STATISTICS,  AND  THEORY  OF  DELAY  CALIBRATION 


DIRECT  ESTIMATION  OF  MEAN  DELAY 
Introduction 

In  radar  and  communication  systems,  it  is  often  desirable  to  know  when  a  pulse  was 
actually  received  at  the  antenna,  not  just  when  it  emerged  from  the  receiver.  Primarily 
because  of  temperature  fluctuations  in  the  receiving  system,  the  delay  through  the  system 
is  constantly  changing.  Thus,  there  is  a  need  for  occasional  delay  calibration  measurements. 

Even  at  constant  temperature,  the  delay  is  a  random  variable,  due  to  the  presence  of 
random  noise  at  the  receiver  input.  Consequently,  the  calibration  procedure  must  estimate 
the  mean  delay.  Usually,  the  time  available  for  calibration  limits  the  number  of  measure¬ 
ments  which  can  be  made  to  estimate  the  mean  delay.  To  determine  the  accuracy  possible 
with  limited  measurements  of  the  delay  is  one  of  the  major  purposes  of  this  report. 

A  method  of  direct  measurement  of  the  delay  is  depicted  in  Fig.  1.  The  technique 
is  to  measure  the  time  interval  between  an  input  mark  pulse  and  the  resulting  output  mark 
pulse.  The  input  pulse  is  generated  by  the  data  simulator;  the  output  pulse  emerges  from 
the  receiver.  The  timing  generator  generates  a  time  of  day  represented  by  a  group  of  bits. 


Fig.  1— Direct  measurement  method  for  delay  calibration 


In  each  timing  distributor,  a  time  of  day  is  assigned  to  the  received  mark  pulse.  The 
difference  between  the  assigned  times  is  computed  digitally  in  the  delay  measurement  unit, 
the  output  of  which  is  the  measured  delay.  Computer  software  uses  multiple  measurements 
to  estimate  the  mean  delay. 

This  method  is  conceptually  simple  and  direct.  It  can  be  highly  automated,  with 
most  of  the  data  recording,  computation,  and  data-presentation  functions  handled  by  the 
computer.  However,  the  very  presence  of  the  measuring  apparatus  disturbs  the  receiving 
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system  and  introduces  inaccuracies  into  the  mean-delay  estimation.  The  primary  dis¬ 
turbance  is  due  to  the  quantizations  introduced  by  the  timing  distributors.  The  effect 
of  the  quantizations  shall  be  thoroughly  analyzed  in  the  next  section. 

Another  source  of  error  is  the  delay  caused  by  the  presence  of  the  modulator  in  the 
calibration  system.  The  modulator,  which  is  only  employed  during  calibration,  generally 
can  be  temperature  controlled  sufficiently  to  insure  that  its  time  of  passage  is  a  constant. 
Thus,  the  modulator  delay  need  only  be  measured  accurately  once.  Subsequently,  this 
constant  can  be  subtracted  from  the  measured  system  delay,  either  by  the  computer 
software  or  in  the  delay  measurement  unit.  Thus,  the  modulator  delay  shall  not  be 
considered  explicitly  in  the  subsequent  sections  of  this  report. 

We  shall  assume  that  the  timing  distributors  are  not  part  of  the  receiving  system  but 
are  inserted  during  calibration  only.  In  each  distributor,  there  is  an  elapsed  time  between 
the  arrival  of  the  pulse  at  the  input  and  the  instant  of  time-of-day  assignment.  If  the 
elapsed  times  for  the  two  distributors  are  not  equal,  a  constant  bias  is  introduced  into 
the  mean-delay  estimation.  The  preceding  discussion  concerning  the  modulator  bias  applies 
to  this  differential  distributor  bias.  Hence,  the  latter  shall  not  be  considered  further. 

In  most  receiving  systems,  the  mean  delay  is  sensitive  to  variations  in  the  signal-to- 
noise  ratio.  Thus,  this  ratio  should  be  held  reasonably  constant  during  calibration.  In  a 
later  section,  a  statistical  test  of  the  sensitivity  shall  be  discussed. 


General  Theorems  for  a  Single  Measurement 

We  shall  define  the  system  delay  d  as  the  sum  of  the  modulator  and  receiver  delays, 
and  the  differential  distributor  delay.  Let  dm  indicate  the  measured  delay,  which  is  defined 
as  the  quantity  actually  measured  and  sent  to  the  computer  or  tape  storage.  As  shown  in 
Fig.  2,  we  can  write 


INPUT  MARK  OUTPUT  MARK 


Fig.  2— Timing  assignments  in  the  direct  measurement  of  delay 


dm  =  d  ~  X1  +  *2 


(1) 


where  Xj  and  x2  are  the  quantization  errors  due  to  the  time  tagging  of  the  input  and 
output  marks,  respectively.  Immediately  following  the  occurence  of  a  mark,  the  time 
tagging  involves  assigning  a  time  of  day  corresponding  to  the  next  leading  edge  of  a  clock 
in  the  timing  generator.  It  is  assumed  that  is  a  random  variable  with  a  uniform  distribu¬ 
tion  and  is  statistically  independent  of  the  random  variable  d.  However,  examination  of 
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Fig.  2  indicates  that  x2  is  completely  specified  by  the  values  of  d  and  jcj  Specifically, 
if  T  is  the  system  clock  period, 


x 


2 


(2) 


where  rem  (z/T)  is  the  remainder  obtained  when  z  is  divided  by  T  if  z  >  0.  For  example, 
if  z  =  16.3T,  then  rem  (z/T)  =  0.3.  When  z  <  0,  we  define  rem  (z/T)  =  rem  [(mT  +  z)/T], 
with  m  such  that  mT  +  z  >  0. 

If  we  define 

y  =  d  -  XV  (3) 


it  follows  that 


=  y  +  T  -  T  rem 

m  ^ 


y 

T 


(4) 


Since  Xj  is  uniformly  distributed  between  0  and  T  and  is  statistically  independent  of 
d,  it  can  be  shown  with  elementary  probability  theory  that  the  density  function  of  the 
random  variable  y  is  given  by 


T 

h(y)  =  -  f  f(y  +  x)dx , 
1  Jo 


(5) 


where  f(z)  is  the  density  function  of  the  random  variable  d.  We  shall  always  assume  that 
f(z)  is  a  well-behaved  continuous  function. 

Observe  that  because  of  the  time-tagging  mechanism,  the  probability  is  zero  that 
dm  is  not  equal  to  a  positive  integral  multiple  of  the  bit  period  T.  F-om  Eq.  (4)  it  is 
seen  that  dm  is  a  function  of  the  single  random  variable  y.  It  follows  from  elementary 
probability  theory  that 


P(dm  =  nT)  =  P(nT  -  T  <  y  <  nT ),  (6) 

where  n  is  a  positive  integer  and  P(  )  denotes  the  probability  of  the  event  in  parentheses. 
Using  Eq.  (5)  in  Eq.  (6),  we  conclude  that  the  probability  of  the  measured  delay  is  given 
by 


P(dm  =  t)  =  0  r  *  nT 

nT  T 

P(dm  =  nT)  =  [  dy  f  dx  ~ f(y  +  x),  n  =  1,  2,  3 .  (7) 

J(n-\)T  J0  1 

If  we  regard  dm  as  an  estimate  of  the  expected  value  of  d,  the  “true”  delay,  it  is 
important  to  establish  whether  the  estimate  is  biased.  Let  E(  )  denote  the  expected  value 
of  the  quantity  in  parentheses.  From  the  definition  of  expected  value  and  Eq.  (7),  we 
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Using  Eq.  (14),  it  follows  that  the  factor  in  brackets  in  Eq.  (13)  is  equal  to  the  constant 
77 2.  From  the  normalization  property  of  density  functions,  the  remaining  integration  if 
trivial.  Collecting  terms,  we  obtain  Eq.  (9),  and  the  theorem  is  proved. 

From  the  definition,  Theorem  1,  and  Eq.  (7),  we  can  express  the  variance  of  the 
measured  delay  in  the  form 


o 


2 

T 


dx  f(y  +  x)  - 


(15) 


This  equation  is  the  moot  convenient  form  for  numerical  computation  once  the  density 
function  is  specified.  However,  for  analytical  purposes  and  to  determine  approximate 
formulas,  an  integral  expression  is  more  useful.  Furthermore,  it  is  desirable  to  obtain 
bounds  on  the  variance  for  an  arbitrary  continue  us  density  function.  The  next  theorem 
provide.,  the  needed  information.  We  denote  the  variance  of  the  system  delay  by  o2. 

Theorem  2.  The  variance  of  the  measured  delay  can  be  expressed  as 


For  any  continuous  density  function,  the  following  inequality  must  be  satisfied: 

o2  <  ol  <  a2  +  — . 

T  4 

Proof:  From  the  expected-value  theorem,  we  write 
E(d2m)  =  |  h(y) 

Using  Eq.  (5)  and  performing  a  change  of  variable,  we  obtain 


y  +  T  -  T  rem  — 
T, 


dy. 


E(d2m)  =  j  f(y)I(y)dy , 

where 

2 

dx. 


Substituting  Eq.  (14)  into  Eq.  (20)  and  rearranging  yields 


i(y) 


-  f 

T 


y  -  x  +  T  -  T  rem 


(16) 


(17) 


(18) 


(19) 


(20) 


5 


/ 


DON  J.  TORRIERI 


i(y) 


-ts: 


y  +  T  -  Trer-i 


Tu  I  x  -  1  rem  |— ) 
T 


dx,  ,21) 


where  u(x)  is  the  unit  step  function  defined  by 

fo,  x  <  0; 

1,  x  >  0. 

L 

Squaring  the  term  in  brackets,  simplifying,  ar  d  noting  that  u2(x)  =  u(x),  we  have 


u(x)  = 


(22) 


I(y) 


y  •  T  -  T  rem 


© 


ir-(- 


r.v> 


T  rem  { — )  |  dx. 


-2 y  -  T  +  2 T  rem 

After  performing  the  remaining  integration  and  regrouping  terms,  we  are  left  with 


(23) 


/(y)  =  y2  +  ^2  rem 


(I)  - 7,2  ren,z  (?) 


(24) 


Since  f(y)  is  the  density  function  for  the  system  delay,  we  know  that 

j  y2f(y)dy  =  o2  +  E2(d). 


(25) 


Combining  Eq.  (19),  (24),  and  (25)  yields 

E(d2J  =  a2  +  E2(d)  +  T2  J  f(y) 


rem 


dy. 


(26) 


From  the  definition  of  o^.,  Theorem  1,  and  Eq.  (26),  we  obtain  the  final  expression  given 
by  Eq.  (16).  The  double  inequality  of  Eq.  (17)  follows  from  Eq.  (16)  by  noting  that 
from  elementary  calculus 


0  <  rem  £)  -  rem2  <  J 


and  employing  the  normalization  property  of  density  functions. 


(27) 


The  variance  o2  is  a  function  of  the  mean  system  delay.  Since  it  can  be  changed  by 
inserting  a  variable  delay  line  into  the  receiving  system,  the  mean  delay  may  be  regarded 
as  an  adjustable  parameter.  We  can  increase  the  accuracy  of  any  estimator  of  the  mean 
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delay  by  setting  the  mean  delay  so  that  a  j,  is  near  its  minimum.  To  this  effect,  we  develop 
two  theorems  which  help  us  to  locate  the  appropriate  mean  delay,  which  we  shall  frequently 
denote  by  y~  =  E(d)  for  simplicity.  Since  we  are  concerned  only  with  physically  realistic 
stituations,  we  restrict  consideration  to  density  functions  of  the  form 

f(y)  =  g(\y  -  yl),  (28) 

where  g(|x|)  is  a  continuously  differentiable  function  of  x. 

Theorem  3.  The  variance  of  the  measured  delay  is  a  periodic  function  of  the  mean 
delay.  The  period  is  equal  to  the  clock  period  T. 

Proof.  Substitution  of  Eq.  (28)  into  Eq.  (16)  and  a  change  of  variables  yields 


o2  =  o2  +  T2 


f  g(ly|)  |rem^-— j  -  rem2(^~^j 


dy. 


(29) 


Clearly,  only  the  factor  in  brackets  is  a  function  of  y.  From  the  definition  of  the  rem 
function, 


rem 


fy  +  y  +  T\  _  (y  +  y\ 

[  T  )  rem\  T  /' 


(30) 


Thus,  the  factor  in  brackets  in  Eq.  (29)  is  periodic  in  y.  It  follows  that  Oj.  is  periodic  in 
y~  and  the  theorem  is  proved. 

Theorem  4.  The  extrema  of  the  variance  of  the  measured  delay  occur  at  the  two  sets 
of  points 


y  =  kT,  k  =  1,2,3,. . .  , 

and 

y  =  kT  +  -,  k  =  1,2,3,...  .  (31) 

St 

The  maxima  occur  at  one  of  these  sets;  the  minima  occur  at  the  other  set. 

Proof.  Under  the  assumption  that  g(bcl)  is  continuously  differentiable,  from  Eqs.  (16) 
and  (28)  we  obtain 

=  T2  r  dg(]y 

37  dy 

Since  g(|y  -y))  is  a  symmetric  function  of  y  -  y",  the  partial  derivative  in  Eq.  (32)  is  an 
antisymmetric  function  of  y  -~y.  At  the  sets  of  points  indicated  in  Eq.  (31),  the  function 
in  brackets  in  Eq.  (32)  is  a  symmetric  function  of  y  -  y.  Consequently,  the  integrand  in 
Eq.  (32)  is  antisymmetric,  and 
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3  o| 

-  =0  (33) 

3y 


at  the  points  defined  by  Eq.  (31).  Thus,  these  points  are  the  extrema.  From  Theorem  3 
and  the  continuity  of  o|,,  one  set  must  be  the  maxima  and  the  other,  the  minima. 

To  determine  which  set  of  points  in  Eq.  (31)  represents  the  minima,  we  can  examine 
the  second  derivative  of  the  variance.  It  is  often  easier  to  evaluate  the  variance  explicitly 
at  each  set  and  determine  the  minimu  by  direct  comparison. 


A  special  case  of  great  interest  occurs  when  the  quantization  is  gross  compared  to 
the  system  standard  deviation  or  when  the  signal-to-noise  ratio  is  high.  In  this  case,  the 
integral  expression  of  Theorem  2  can  be  simplified.  We  assume 

f(y)  >  0,  E(d)  -  T  <  y  <  E(d)  +  T\ 

l\y)  585  0,  otherwise.  (34) 

We  write  E(d)  in  the  form 


E(d)  =  KQT  +  R 


(35) 


where  KQ  is  the  integral  part  obtained  when  E{d)  is  divided  by  T,  and  R  is  defined  by 


R  =  T  rem 


(36) 


If  we  use  Eqs.  (34)  and  (35),  Theorem  2  yields 

C  K0T+R+T 

o\  =  O2  +  T2  |  f(y ) 

JK0T+R-T 

The  integral  can  be  divided  into  three  regions.  In  these  regions  we  have 


rem 


(|)  -  (|) 


dy. 


(37) 


fy  -  V 

T 


I  y  -  k0t 
l  T 

On  substitution  of  Eq.  (38)  into 


+  1,  K0T-T  +  R<y<K0T\ 

K0T<y<K0T  +  T-, 

-  1,  K0T  +  T<y<K0T+T  +  R. 

Eq.  (37)  and  some  manipulation,  there  results 


(38) 
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77 


it  ••  a2  +  T2 


+  T 


L 


K0T+R+T 


f(y) 


K.0T*R-T 


y  -  K0T 


dy 


rK  oT  f-2  y  +  2  KnT\ 

j  fiy)[  ■  -f  — ) 

^  k'  Ti.D_T  \  J  / 

T2  f  °  f(y) 

•'IT  Ti.T  \ 


fly)  -2  + 


dy 

2y  -  2KqT 


k0t+t 


dy. 


(39) 


In  view  of  Eq.  (34),  we  can  extend  the  limits  of  the  first  integral  in  Eq.  (39)  to 
infinity.  Making  use  of  Eqs.  (12)  and  (25)  and  the  normalization  condition,  the  first 
two  terms  in  Eq.  (39)  can  be  simplified.  From  Eq.  (28)  and  a  change  of  variables,  the 
last  two  integrals  in  Eq.  (39)  can  also  be  simplified.  Our  final  result  is 


CT  rT 

at  =  R(T  -  R)  +  9.T  (y  -  R)g(\y\)dy  +  2 T  (y  ••  T  +  «)g(lyl)dy. 

JR  JT-R  (40) 

Note  that  is  symmetric  about  R  =  77 2.  From  the  above,  it  is  apparent  that 

-T 


o2(R  =  0)  =  2 T  f  yg(\y\)dy. 
J(\ 


(41) 


According  to  Theorems  3  and  4,  Eq.  (41)  gives  the  value  of  a2  at  one  set  of  extrema.  In 
view  of  Eqs.  (28)  and  (34),  we  can  write 

ot(kT)*TE[ ly  -  yll,  k  =  1,  2, 3 .  (42) 

As  we  shall  see  in  the  next  section,  this  expression  usually  gives  the  minimum  variance. 
Evaluating  Eq.  (40)  at  R  =  T/2  and  using  the  inequality  in  Theorem  2,  we  obtain  the 
bounds: 

T2  /  T\  T2 

—  <  4  \kT  +  +  °2,  *  =  1,2,3 .  (43) 


This  variance  is  generally  the  maximum. 


Normally  Distributed  System  Delay 

In  this  section,  we  obiain  formulas  for  the  probability  distribution  and  the  variance 
under  the  assumption  that  the  system  delay  is  approximately  described  by  a  normal 
density  function.  Thus,  we  take 


1 

f(y)  =  -  exp 

y/T.iO 


(y  -  y)2 

2o2  . 


(44) 
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where  ~y  =  E(d)  and  o2  is  the  variance  of  the  system  delay.  Since  the  system  delay  can 
never  take  a  negative  value,  we  must  assume  "y  »  o  to  use  Eq.  (44)  for  all  y.  Substituting 
into  Eq.  (7),  we  have 


where 


1 

erfc  (jc)  = - 


dx 


(45) 


(46) 


is  the  complementary  error  function.  By  a  suitable  change  of  variable  in  each  term,  Eq. 
(45)  can  be  written  in  the  form, 


P(dm  =  nT)  =  ^  f  eric  y  dy  - 
1  Ja-Tlo 


-  f 

T  J. 


a*Tlo 


erfc  y  dy. 


where  for  convenience  we  have  defined 

nT  -  y 

a  =  - . 

a 

From  tabulated  integrals,  it  can  be  shown  that 

ferfcxdhc  =  xerfex  -  N(x), 

where  for  convenience  we  have  defined  the  standard  normal  function, 


(47) 


(48) 


(49) 


N(x)  = 


exp 


(50) 


Using  Eq.  (49)  in  Eq.  (47),  we  obtain  a  closed-form  expression  for  the  probability  distribu¬ 
tion  of  the  measured  delay  when  the  system  delay  has  a  normal  density  function: 


P(dm  -  r)  -  0, 


r  ¥=  nT; 
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In  tvro  special  cases,  approximate  closed-form  expressions  can  be  obtained  for  the 
variance  of  the  measured  delay.  First  consider  the  case  where  the  quantization  is  fine. 
Thus,  the  system  standard  deviation  is  much  greater  than  the  clock  period,  that  is,  a  »  T 
We  employ  the  Taylor  series  expansions  valid  for  small  values  of  b: 

Ni't  +  b)  *  JV(a)[l  -  aft]; 

erfc  (a  +  b)  *  erfca  -  bN(a).  (52) 

Substituting  Eq.  (52)  into  Eq.  (51)  yields  the  approximate  relation 

P(dm  =  nT)  =  -N(a),  o  »  T,  n  -  1,2,3 .  (53) 

m  a 

Thus,  in  this  case  the  limiting  form  of  the  envelope  of  the  probability  distribution  is 
normal.  It  follows  that 

otp  -*  o2,  a  »  T,  (54) 

which  is  intuitively  satisfying. 

The  second  special  case  to  be  considered  occurs  when  Eq.  (34)  is  applicable.  For  a 
normal  distribution,  the  condition  T  >  4a  will  be  accepted  as  sufficient  for  Eq.  (34).  We 
can  then  employ  Eq.  (40),  with  g(y)  =  N(y).  From  a  tabulated  integral,  it  can  be  shown 
that 


jxN(x)dx  =  -N(x)  (55) 

Using  this  equation  and  omitting  terms  which  are  negligible  for  all  values  of  R,  we  obtain 

o \  *  R(T  -  if)  -  2 T(T  -  if)  erfc  ~  ^  -  2TR  erfcf^j 


+  2  To 


'<-¥)  •  “(4) 


T  >  4o. 


(56) 


Taking  the  second  derivative  of  this  expression  with  respect  to  if,  we  get 

8 2o2t  _  2 T 
8if 2  ~  ~a 


(T  -  R\  (R\ 

N\— j*N[-j 


-  2. 


(57) 


By  direct  substitution,  it  follows  that  a|,  is  a  minimum  when  if  =  0  and  a  maximum 
when  if  =  T/2.  Omitting  negligible  terms,  Eq.  (56)  and  Theorem  4  imply 


Minima:  o^(feT) 


To,  k  =  1,  2, 3,  . . .  ; 


(58) 
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Maxima:  o|. 


kT  +  *  -L  +  \ToN  -  2 T2  erfc  T  >  4a  , 


k  =  1,  2,  3 . 


(59) 


Of  course,  Eq.  (58)  could  have  been  derived  directly  from  Eq.  (42). 


Multiple  Measurements 

During  the  period  of  calibration,  many  delay  measurements  are  made.  If  the  mein 
delay  is  regarded  as  an  unknown  but  nonrandom  parameter,  one  might  wish  to  make  a 
maximum-likelihood  estimate  based  on  the  sa.nple  data.  However,  even  if  the  appropriate 
probability  density  function  is  known,  the  amount  of  computation  required  for  the 
maximum-likelihood  estimate  is  prohibitive.  For  example,  if  the  system  delay  is  normally 
distributed,  the  measured  delay  has  a  single  sample  density  function  given  by  Eq.  (51). 
Gearly,  the  differentiation  of  the  corresponding  multiple  sample  density  function  and  the 
solution  of  the  likelihood  equation  is  a  formidable  task. 

Under  the  circumstances,  a  logical  and  simple  procedure  is  to  use  the  sample  mean 
to  estimate  the  mean  delay.  If  we  let  the  symbol  dm  denote  the  sample  mean  and  dm  ■ 
denote  a  sample  value  of  the  measured  delay,  we  have 

^  ,60> 

i-1 


Looking  at  Fig.  2,  it  is  seen  that  the  samples  can  be  assumed  to  be  independent  if  the 
interval  between  input  marks  is  much  larger  than  the  system  clock  period  T  and  if  the 
system  clock  is  independent  of  the  data  simulator  clock.  If  the  temperature  and  the  signal- 
to- noise  ratio  are  approximately  constant  during  calibration,  we  can  accurately  assume  that 
E(dmj)  =  E(dm ).  Under  thes'  conditions,  it  follows  from  Eq.  (60)  and  Theorem  1  that 

E(dJ  =  E(dm)  =  E(d).  (61) 


Similarly,  we  can  write  the  following  expression  for  the  variance  of  the  sample  mean: 


a 


2 

m 


(62) 


If  the  number  of  sample  data  points  N  is  large,  we  can  invoke  the  central  limit  theorem 
and  assert  that  the  sample  mean  is  normally  distributed. 

Equation  (62)  indicates  that  we  should  make  a  £  as  small  as  possible.  With  a  variable 
delay  line,  we  can  adjust  the  mean  delay  so  that  o|,  is  a  minimum.  If  the  system  delay  is 
normally  distributed,  it  follows  from  Eqs.  (58)  and  (62)  that  we  can  adjust  the  mean  delay 
so  that 


T  >  4o  . 


(63) 


NRL  REPORT  7812 


As  an  example,  suppose  T  =  100  ns,  a  =  20  ns,  N  =  1000,  and  d  =  10.173  /is.  With 
the  insertion  of  a  delay  of  27  ns  into  the  receiving  system,  we  obtain  dm  =  10.2  /is,  which 
is  an  integral  multiple  of  T.  Our  variance  is  now  approximately  a  minimum,  and  from 
Eq.  (63)  we  know  that  am  ^  1.3  ns. 


If  the  probability  distribution  of  the  system  delay  is  unknown  or  possesses  unknown 
parameters,  the  variance  of  the  measured  delay  cannot  be  calculated  using  Eq.  (42).  We  shall 
discuss  the  empirical  estimation  of  the  probability  distribution  in  the  next  section.  Another 
approach  to  determining  o £  is  to  use  the  estimate  provided  by  the  sample  variance,  which  is 
defined  by 

S2  ■  iTTlIX.  (64) 

(“I 

Under  the  same  assumptions  leading  to  Eq.  (61),  we  have 

E(S2)  =  a (65) 

Thus,  the  sample  variance  provides  an  unbiased  estimate  of  o|,.  It  follows  from  Eq.  (62) 
that  an  unbiased  estimate  of  the  variance  of  the  sample  mean  is  given  by 


m 


S2 
~N  ' 


(66) 


where  the  circumflex  above  the  symbol  on  the  left-hand  side  indicates  an  estimated 
quantity. 

To  precisely  determine  the  accuracy  provided  by  the  estimate  of  Eq.  (66),  we  need  to 
know  the  probability  distribution  of  the  measured  delay.  However,  an  approximate  computa¬ 
tion  of  the  variance  of  the  estimate  can  be  accomplished  by  assuming  that  dm  is  approximately 
normally  distributed.  It  then  follows  that  ( N  -  1  )S2/Oy  has  a  chi-square  distribution  with 
N  -  1  degrees  of  freedom.  From  this  fact  and  Eq.  (66),  we  obtain  the  variance  of  the 
estimated  variance: 


var  (o2  ) 

'  m ' 


20* 

N^(N  -  1)' 


(67) 


Another  use  for  the  sample  variance  is  in  estimating  the  standard  deviation  of  the 
system  delay,  which  can  often  be  related  theoretically  to  the  signal-to-noise  ratio  at  the 
receiver  input.  If  the  system  delay  has  a  normal  distribution  and  the  mean  delay  is 
appropriately  adjusted,  Eq.  (58)  is  valid.  If  the  left-hand  side  of  Eq.  (58)  is  estimated  by 
the  sample  variance,  we  obtain  the  following  estimate  of  the  standard  deviation  of  the 
system  delay: 


a  = 


IF' 


The  same  crude  approximation  used  to  obtain  Eq.  (67)  now  yields 


(68) 
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var  (o) 


no%, 

T^IN  -  1) ' 


(69) 


It  was  mentioned  in  the  first  section  that  the  mean  delay  is  sensitive  to  variations  in 
the  signal-to-noise  ratio.  Thus,  we  would  like  to  know  over  what  range  of  signal-to-noise 
ratio  the  mean  delay  can  be  considered  a  constant.  Suppose  that  we  measure  the  sample 
means  at  various  fixed  values  of  signal-to-noise  ratio.  Let  aTj  and  x2  be  the  sample  means 
from  two  random  samples  of  size  n1  and  n2  at  two  different  signal-to-noise  levels. 
Assuming  that  the  two  random  samples  are  independent,  we  test  the  null  hypothesis  that 
the  two  population  means  are  equal  against  the  alternative  that  they  are  not.  If  n{  and 
n2  are  large  enough  that  we  can  invoke  the  central  limit  theorem,  the  test  statistic 


Z  = 


(70) 


has  a  standard  normal  distribution,  where  Oj  and  o2  are  the  population  variances.  To 
test  the  null  hypothesis  at  the  confidence  level  a  =  0.05,  the  critical  region  is  Z  >  1.96 
When  the  variances  are  unknown,  we  substitute  the  sample  variances  S ^  and  foroj 
and  o2.  We  then  accept  the  null  hypothesis,  if 


I*, 


-  x2  j  <  1.96 


(71) 


Repeating  this  test  over  all  pairs  in  a  range  of  signal-to-noise  ratio,  we  can  validate  or 
reject  the  approximate  constancy  of  the  mean  delay  within  the  range. 


EMPIRICAL  DISTRIBUTIONS 
Manual  Distribution  Method 

To  make  use  of  the  formulas  and  theorems  derived  in  the  prior  sections,  we  must  have 
at  least  some  knowledge  of  the  density  function  of  the  system  delay.  For  example. 
Theorem  4  can  be  invoked  only  if  it  can  be  verified  that  f(y)  is  symmetric  about  the 
mean,  that  is,  if  Eq.  (28)  can  be  shown  to  hold.  We  shall  now  examine  the  empirical 
estimation  of  the  probability  distribution  function. 

The  manual  distribution  method  (MDM)  is  implemented  as  shown  in  Fig.  3.  It  is 
assumed  that  the  input  marks  constitute  a  stationary  ergodic  process  during  the  time 
period  of  measurement.  The  two  counters  are  set  to  count  for  a  fixed  time  interval. 

The  gate  pulse  of  the  digital  time  delay  unit  is  triggered  by  the  occurrence  of  the  input 
mark.  Throughout  the  duration  of  the  gate  pulse,  the  high-speed  switch  will  not  pass 
any  output  marks.  Any  output  marks  received  after  the  fall  time  of  the  gate  pulse 
will  pass  through  the  high-sp  ;ed  switch.  Initially  the  duration  of  the  gate  pul  e,  r, 
is  set  to  the  smallest  value  at  which  the  number  of  events  registered  on  the  passed-pulse 
counter,  N  ,  is  less  then  the  total  number  of  marks,  N,  registered  on  the  mark 
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Pig.  3— Instrumentation  of  the  manual  distribution  method 


counter.  The  setting  r  is  then  stepped  in  short  intervals  for  successive  measurements  of 
N^.  Each  measurement  of  Np  determines  a  point  of  the  empirical  distribution  function, 
given  by 


F(r)  =  1  -  ^  ■  (72) 

In  other  words,  F(t)  is  Em  estimate  of  the  probability  that  the  measured  delay  is  less  than 
r.  Linear  extrapolation  ordinarily  is  used  to  determine  F(t)  for  values  of  r  between  settings. 

The  MDM  is  a  time-consuming  procedure.  To  insure  ergodicity,  temperature  fluc¬ 
tuations  and  the  signal-to-noise  ratio  at  the  receiving  system  input  must  be  carefully  con¬ 
trolled  while  the  MDM  is  being  employed.  Furthermore,  the  time  available  for  calibration 
im  /  be  insufficient  for  application  of  the  MDM  to  an  operating  system.  For  these  reasons, 
the  technique  is  useful  primarily  during  design,  initial  deployment,  Euid  malfunction 
correction. 

We  shall  define  the  system  delay  as  the  delay  due  to  the  modulator  and  receiver. 

The  quantity  which  actually  is  compared  to  the  setting  r  shall  be  called  the  measured  delay. 
The  MDM  produces  an  empirical  distribut  on  function  for  the  measured  delay.  ,An  estimate 
of  the  mean  value  of  the  measured  delay  is  provided  by  the  value  rQ  at  which  F(tq)  =  1/2. 
From  the  empirical  distribution  function,  the  corresponding  density  function  can  be 
calculated  numerically.  A  numerical  integration  csin  then  be  performed  to  derive  various 
moments,  such  eis,  the  variance,  skewness,  and  kurtosis.  We  shall  examine  next  the 
accuracy  of  the  MDM  and  establish  a  model  for  the  relation  of  the  measured  delay  to  the 
system  delay. 

Because  the  high-speed  switch  cannot  respond  instantaneously  to  the  gate  pulse, 
there  is  a  time  lag  between  the  fall  time  of  the  gate  pulse  and  the  actual  time  at  which 
the  output  mark  pulses  begin  to  pass  through  the  switch.  Consequently,  a  constant  bias 
is  introduced  into  the  measured  delay.  The  rise  time  of  the  gate  pulse  coincides  with  a 
clock  pulse  internal  to  the  digital  delay  unit.  We  shall  denote  this  clock  period  by  T. 

If  the  clock  in  the  data  simulator  is  nonsynchronous  with  the  digital  delay  unit  clock  and 
if  the  marks  are  separated  by  many  clock  periods,  the  uncertainty  in  the  start  time  of  the 
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gate  pulse  can  be  considered  a  uniformly  distributed  quantization  error,  which  ranges  from 
0  to  T  seconds.  Thus,  we  have  established  the  following  model: 

dm  =  d  -  c  +  *  *  <73> 

where  dm  is  the  measured  delay,  d  is  the  system  delay,  c  is  the  bias  of  the  high-speed 
switch,  and  x  is  the  quantization  error. 


Mathematical  and  Experimental  Analysis  of  the  MDM 

Since  the  random  variables  d  and  x  are  independent  and  x  is  uniformly  distributed, 
it  follows  from  Eq.  (73)  that 


and 


E(d„)  =  E(d)  -  c 


(74) 


(75) 


where  o2  is  the  variance  of  dm  and  o2  is  the  variance  of  d.  Equation  (74)  tells  us  how  to 
compute  the  mean  system  delay,  which  is  what  we  are  really  trying  to  estimate,  from  the 
mean  measured  delay,  which  is  estimated  from  the  empirical  distribution  function  supplied 
by  the  MDM.  Equation  (75)  indicates  that,  generally  speaking,  the  shape  of  F(t)  will  not 
be  significantly  affected  by  the  quantization  error  if  a  is  at  least  several  times  the  value  of 
T. 


A  test  of  our  model  and  an  evaluation  of  the  bias  can  be  accomplished  by  the 
experimental  procedure  illustrated  in  Fig.  4.  The  two  pulse  generators  are  activated  by 
the  same  external  clock.  Thus,  the  simulated  input  and  simulated  output  pulses  have  a 
fixed  difference  in  starting  times,  which  we  shall  call  the  simulated  delay.  Because  of  the 
common  external  clock,  the  simulated  delay  must  be  an  integral  multiple  of  the  clock 
period.  If  the  estimated  mean  value  of  the  delay  is  not  an  integral  multiple  of  the  clock 
period,  the  deviation  provides  a  measure  of  the  bias.  If  the  external  clock  is  stable  during 


Fig.  4— Model  test  and  bias  measurement 
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the  test,  the  simulated  delay  will  have  negligible  variance.  Thus,  F(t)  will  show  the 
distribution  of  the  quantization  error. 

As  an  example,  Fig.  5  shows  some  actual  data  taken  from  a  model  test  using  a 
Tektronics  7D11  Digital  Time  Delay  Unit.  This  unit  has  an  internal  clock  with  a  2-ns 
period.  The  simulated  input  pulses  were  separated  by  200  ps,  and  we  let  N  =  34,247. 

In  Fig.  5  the  straight  line  indicates  a  uniform  distribution  over  a  2-ns  interval.  The  actual 
data  is  seen  to  be  quite  close  to  the  ideal. 

The  external  clock  feeding  the  pulse  generators  had  a  period  of  1  ps,  and  the  simulated 
delay  was  45  ps.  As  shown  in  the  figure,  the  mean  measured  delay  is  estimated  as  44.9343 
ps.  Using  Eq.  (74),  we  arrive  at  the  estimated  value  of  the  bias  as  c  =  66.7  ns. 


TIME  DELAY  (mcrosecondj) 

Fig.  5— Data  for  the  MDM  model  test 


For  any  fixed  value  of  t,  F(t)  is  a  random  variable  which  can  take  only  rational  values. 
From  Eq.  (72),  the  probability  of  any  particular  rational  number  is 


F(t)  =  — 
N 


=  P[Np  =  N 


k],  k  =  0,  1,  2,  .  .  .  N, 


(76) 


where  P  (  ]  denotes  the  probability  of  the  event  enclosed  in  the  brackets.  The  number 
Np  is  determined  by  successive  Bernoulli  trials  of  the  measured  delay.  Hence,  N  has  a 
binominal  distribution.  Letting  1'(t)  denote  the  “true”  distribution  function  of  the  measured 
delay,  Eq.  (76)  then  yields 
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nr)  =  - 


From  this  relation,  we  obtain 


■0 


lF(r)]k[  1  -  F(r)]^-* 


(77) 


£[F(r))  =  F(t) 


(78) 


and 


2  F(r)[l  -  F(r)] 

ot(r)  =  - 

F  N 


(79) 


where  the  leftrhand  side  of  Eq.  (79)  is  the  variance  of  F(t)  at  a  fixed  value  of  t. 
Equation  (78)  states  that  the  MDM  provides  an  unbiased  estimate  of  F(t)  at  each  value 
r.  Equation  (79)  gives  a  measure  of  the  deviation  of  F(r)  with  respect  to  F(t).  Since 
0  <  F(r)  <  1,  we  obtain  the  upper  bound  given  by 


9  1 

a?(r)  <  — , 
Fx  4 N 


(80) 


where  equality  holds  at  the  value  Tj  for  which  F(Tj)  =  1/2. 

We  now  examine  the  accuracy  with  which  the  mean  measured  delay  can  be  estimated 
from  the  empirical  distribution  function.  A  rigorous  derivation  involves  the  adaptation  of 
arrival-time  estimation  theory  [3].  The  derivation  is  so  complicated  that  it  deserves  a 
report  by  itself.  Hence,  we  shall  content  ourselves  with  an  approximate  calculation  of  the 
variance  of  the  mean  delay  estimated  from  F(r). 

We  assume  that  the  resolution  and  the  number  of  settings  of  r  are  such  that,  with  negli¬ 
gible  error,  F(r0)  =  1/2  for  some  setting  t0  of  the  digital  time  delay  unit.  In  other  words, 
we  exclude  the  possibility  that  F(tq)  =  1/2  for  an  extrapolated  setting.  If  the  number  of 
points  N  is  large,  then  it  follows  that,  for  the  setting  rQ,  we  have  essentially  measured 
the  median  of  the  sample.  We  can  then  apply  the  following  theorem  [1] ,  which  holds 
when  the  population  density  fly)  is  continuous  and  nonzero  at  the  population  median  p: 

For  large  n,  the  sampling  distribution  of  the  median  for  random  samples  of  size  2n  +  1 
is  approximately  normal,  with  the  mean  p  and  the  variance  l/8n[f(j£f)] 2.  If  f{y)  is 
symmetric  about  the  population  mean  p  =  E(dm ),  then  P~P-  Assuming  that  N  is  sufficiently 
large  that  N  =  2n,  we  conclude  that  the  estimated  mean  measured  delay  is  approximately 
normal,  with  variance 


4N[fl/i)]2' 


(81) 


If  f(p)  is  unknown,  we  can  estimate  it  by  numerical  differentiation  of  F(p).  We  can  test 
whether  two  estimated  means  are  equal  by  using  Eq.  (71),  with  the  estimated  mean  in 
place  of  x  and  o?  in  place  of  S2/n. 
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Nonlinear  Regression  Analysis 

When  the  “true”  distribution  function  of  the  measured  delay  is  known,  except 
for  some  parameters,  these  parameters  can  be  estimated  from  the  empirical  distribution 
function.  If  the  parameters  are  expressed  in  terms  of  the  moments  of  the  density  function, 
they  can  be  estimated  by  the  numerical  procedure  described  in  the  MDM  section.  Alterna¬ 
tively,  the  pa  ameters  can  be  estimated  by  the  least-squares  method  of  nonlinear  regression 
analysis.  This  method  entails  the  estimation  of  the  parameters  by  minimizing  the  sum  of 
squares  of  deviations  of  F(t,)  with  respect  to  F(Tt),  where  the  r(  are  the  values  of  the 
discrete  settings. 

A  more  common  situation  is  when  one  wishes  to  test  a  proposed  regression  model. 

The  parameters  of  the  model  can  be  determined  by  the  numerical  procedure  or  the  least- 
squares  method.  Then,  the  fit  of  the  model  is  tested  by  an  appropriate  “goodness  of  fit” 
statistic. 

The  standard  goodness  of  fit  tests,  such  as  the  chi-square  and  the  Kolmogorov- 
Smimov  [2],  require  a  data  format  not  supplied  by  the  MDM.  Thus  special  tests  must  be 
applied.  We  shall  investigate  two  different  techniques,  both  of  which  use  the  statistics 


F( Tj )  -  F(rf) 

Z:  =  VAT  -  ,  i  «  1,  2,  . .  .  n,  (82) 

where  n  is  the  number  of  discrete  settings  r(.  which  are  such  that  F(rf)  0,  1.  As  defined 
previously,  N  is  the  number  of  samples  per  fixed  setting.  It  shall  be  assumed  that  F(r,)  is 
independent  of  F(t.)  when  j  #  i.  Thus,  the  Z(are  independent  of  each  other.  If  N  is 
sufficiently  large,  the  binomial  distribution  of  F(r-)  can  be  approximated  by  a  normal  one. 
From  Eqs.  (78)  and  (79),  it  then  follows  that  the  Z,  have  standard  normal  distributions. 
The  first  test  of  the  null  hypothesis  that  F(r)  is  the  true  distribution  function  is  an  im¬ 
mediate  consequence  of  these  considerations. 

Modified  Chi-Square  Test: 

Form  the  test  statistic 


X2  =  £  Zf  (83) 

r=l 


Test  this  statistic  as  a  chi-square  distribution  with  n  -  q  degrees  of  freedom,  where  q  is 
the  number  of  parameters  estimated  from  F(t). 


The  second  test  is  based  on  the  simple  fact  that  either  all  of  the  |Z(|  are  less  than  a 
fixed  number  D  or  at  least  one  of  the  |Z.|  is  greater  than  D.  Thus,  we  may  write 


P{ some  |Z.|  >  D)  =  1 
=  1 
=  1 


P(all  |Z(.|  <  D) 
[P(IZl<D)]" 

(1  -  2  erfc  D)n , 


(84) 
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where  the  second  equality  follows  from  the  independence  of  the  Zi  and  the  last  equality 
results  from  using  Eqs.  (44)  and  (46).  Equation  (84)  provides  the  justification  for  the 
following  test. 

Modified  Kolmogorov-Smimov  Test: 

Determine  the  maximum  absolute  value  of  the  Zt\  that  is,  determine  the  quantity 

D  =  max  |Zf| .  (85) 

i 

Compute  the  test  value 


p  =  1  -  (1  -  2  erfc  D)n.  (86) 

The  null  hypothesis  is  acceptable  at  the  confidence  level  a  if  p  >  a.  Otherwise,  the  null 
hypothesis  is  rejected. 

As  mentioned  previously,  we  can  usually  assume  that  the  measured  delay  is  distributed 
the  same  as  the  system  delay  if  o  is  at  least  several  times  the  value  of  T.  If  the  latter  is 
not  the  case,  we  must  calculate  the  known  or  proposed  distribution  of  the  measured  delay 
by  means  of  Eq.  (73).  If  f(y)  is  the  density  function  of  d  -  c  and  x  is  uniformly  distributed 
from  0  to  T,  it  follows  from  elementary  probability  theory  that 


If  the  system  delay  is  normally  distributed,  we  can  substitute  Eq.  (44)  into  Eq.  (87) 
to  obtain 


F(t)  =  1  -  -  [a  erfc  a  -  N{a ) 


(a  -  b)  erfc  (a  -  b)  +  N(a  -  6)] , 


(88) 


where  we  define 


a 


t  -  E(d)  +  c 
a 


and 


T 

b  =  -.  (89) 

o 


Since  F(p)  =  1/2  when  p  =  E(dm )  =  E(d)  -  c  +  T/2,  the  median  of  the  measured  delay  is  equal 
to  its  mean.  Differentiating  Eq.  (88)  with  respect  to  t,  we  obtain  the  density  function, 


1 

f(r)  =  —[erfc  (a  -  b)  -  erfc  a].  (90) 


20 


NHL  REPORT  781 2 


From  this  relation  and  Eq.  (81),  we  have 

,  2  - 


y2 


0‘ 


4/V 


1-2  erfc  — 
2o 


for  the  variance  of  the  numerically  estimated  mean  measured  delay. 


(91) 


Recent  Developments 

With  the  production  of  lewlett  Packard’s  HP  5345A  reciprocal  counter  [4],  highly 
precise  delay  calibration  is  possible  with  a  single  measuring  instrument.  Whereas  the 
conventional  counter  measures  the  number  of  input  events  during  an  interval,  the  reciprocal 
counter  measures  the  time  interval  between  events. 

The  reciprocal  counter  has  two  basic  modes  of  operation.  In  the  time-interval- 
averaging  mode,  it  may  be  used  in  place  of  the  timing  distributors,  the  timing  generator, 
and  the  delay  measurement  unit  of  Fig.  1.  The  sample  mean  and  sample  variance  are 
determined  according  to  Eqs.  (60)  and  (64).  A  single  delay  measurement  can  be  modeled 
by  Eq.  (1).  However,  band-limited  noise  is  added  to  the  counter’s  time  base,  causing  Xj 
and  x 2  to  be  independent  and  uniformly  distributed  over  the  counter’s  clock  period,  which 
is  T  =  2  ns.  The  use  of  the  noise-modulated  clock  removes  any  harmonic  relationship 
between  the  counter’s  clock  frequency  and  the  repetition  rate  of  the  input  pulses,  helping 
to  insure  the  independence  of  successive  samples  of  ihe  measured  delay.  Noting  that  we 
now  have 


(92) 


Eqs.  (61),  (62),  (65),  (66),  and  (67)  are  valid. 

In  another  mode  of  operation,  depicted  in  Fig.  6,  the  counter  is  used  in  conjunction 
with  a  calculator  and  a  plotter.  The  sample  mean  and  sample  variance  are  determined 


Fig.  6  — Use  of  the  reciprocal  counter 
in  delay  measurement 
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according  to  Eqs.  (60)  and  (64).  A  single  delay  measurement  can  be  modeled  exactly  as 
in  the  direct  method  of  the  first  part  of  this  report.  Hence,  Theorem  2  gives  bounds  on 
the  variance  of  the  measured  delay.  This  mode  is  an  alternative  to  the  manual  distribution 
method.  Each  sample  value  of  the  measured  delay  is  sent  to  the  calculator,  where  it  is 
assigned  to  one  of  a  set  of  disjoint  class  intervals.  In  view  of  Eq.  (7),  the  partition 
bet  veen  adjacent  intervals  should  be  located  halfway  between  the  possible  values  of  the 
measured  delay,  which  are  integral  multiples  of  T  =  2  ns.  After  all  the  sample  values  are 
assigned  to  intervals,  the  calculator  and  plotter  assemble  a  histogram  (empirical  probability 
distribution).  The  conventional  chi-square  or  Kolmogorov-Smimov  tests  can  be  applied 
to  this  histogram  to  determine  the  goodness  of  fit  of  a  proposed  probability  distribution. 
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